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We study how varying the Peclet number (Pe) affects the steady state sedimentation of colloidal 
particles that interact through short-ranged attractions. By employing a hybrid molecular dynamics 
simulation method we demonstrate that the average sedimentation velocity changes from a non- 
monotonic dependence on packing fraction <j> at low Pe numbers, to a monotonic decrease with <j> at 
higher Pe numbers. At low Pe number the pair correlation functions are close to their equilibrium 
values, but as the Pe number increases, important deviations from equilibrium forms are observed. 
Although the attractive forces we employ are not strong enough to form permanent clusters, they 
do induce transient clusters whose behaviour is also affected by Pe number. In particular, clusters 
are more likely to fragment and less likely to aggregate at larger Pe numbers, and the probability 
of finding larger clusters decreases with increasing Pe number. Interestingly, the life-time of the 
clusters is more or less independent of Pe number in the range we study. Instead, the change in 
cluster distribution occurs because larger clusters are less likely to form with increasing Pe number. 
These results illustrate some of the subtleties that occur in the crossover from equilibrium like to 
purely non-equilibrium behaviour as the balance between convective and thermal forces changes. 
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I. INTRODUCTION 

The Peclet number (Pe) measures the ratio of convec- 
tive to thermal forces. At low Pe numbers, where thermal 
forces dominate, one might expect that properties such 
as the structure of a colloidal suspension would be well 
approximated by their values at equilibrium (where Pe is 
strictly zero). On the other hand, at larger Pe numbers, 
purely non-equilibrium phenomena should become more 
prominent [l|, Q ■ 

Nonetheless, the precise conditions under which 
the transition from equilibrium like to purely non- 
equilibrium behaviour occurs is a subtle question and de- 
pends on which properties are under investigation. For 
example, equilibrium correlation functions are often ap- 
plied to calculate the properties of non-equilibrium sys- 
tems. At what Pe number does this approximation break 
down? 

To address such questions, we study the interplay of 
thermal and convective forces in steady-state sedimenta- 
tion, a classical non-equilibrium problem that, despite its 
apparent simplicity, it is still far from completely under- 
stood [ll-Q- A major difficulty for theories arises from 
the long-ranged (1/r) nature of the solvent induced hy- 
drodynamic interactions (HI) that couple the motion of 



the particles in a complex way. Indeed, although Stokes 
[H calculated the sedimentation velocity v® of a single 
hard sphere in 1851, it was necessary to wait more than 
120 years until Batchelor [f| derived the first correction 
in the dilute limit, which is given by v s /v® = 1 — 6.55(f), 
where v s is the average sedimentation velocity of the col- 
loids and <f> is the particle volume fraction. Although this 
is strictly a non-equilibrium effect, it has been shown by 
simulations Q that for hard sphere (HS) particles, the 
behaviour of v s (<p) is virtually independent of Pe num- 
ber down to at least Pe — 0.1. So, for this example, a 
non-equilibrium effect persists virtually unchanged well 
into a regime where the thermal forces are significantly 
stronger than the convective forces. 

In this paper, we introduce a thermodynamic compo- 
nent to this non-equilibrium problem by studying the 
effect of short-ranged attractive interactions on colloidal 
sedimentation at different Pe numbers. Nearly 30 years 
ago Batchelor calculated the effect of short-range interac- 
tions beyond the HS model on the average sedimentation 
velocity in the dilute limit 0, finding: 

v s /v° s = 1 - [6.55 - 3.52(1 - B%)] 4> + 0{<p 2 ) (1) 



where B 2 



Bo is the second virial coefficient, 
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and B2 is the virial coefficient calculated with the ef- 
fective HS radius of the colloids. Eq. (pj suggests that 
adding attractions should increase the sedimentation ve- 
locity, while adding repulsions should decrease it com- 
pared to the pure HS case. An intuition for the increase 
in sedimentation velocity with attractions at low pack- 
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ing fractions can be obtained from the following consid- 
erations: Attractions should increase the probability of 
(transient) cluster formation. The sedimentation force 
scales linearly with the number i of particles in the clus- 
ter, whereas the friction is proportional to the radius 
of gyration of the cluster, which typically scales as i a , 
so that the sedimentation velocity of a cluster scales as 
v s ~ i 1_Q . For compact clusters a rs | whereas for 
more open clusters a can be larger, although it cannot 
be greater than one. Even in the regime where the attrac- 
tions are not strong enough to form permanent clusters, 
transient clusters can form, leading to a speed up in the 
average sedimentation velocity. 

In the dilute limit, the increase in v s with increas- 
ing strength of attractions has been observed experimen- 
tally d, Q. On the other hand, as the packing frac- 
tion increases, the gap between nearest neighbours de- 
creases and so the fluid locally experiences an enhanced 
friction, which eventually slows down the average sedi- 
mentation velocity. What happens in between these two 
limiting regimes has been recently studied by Gilleland 
et al. [10|, ELL! , who calculated v s (0) theoretically in the 
Pe C 1 limit, and who first predicted a maximum as a 
function of </> for strong enough attractions. We applied 
the same mesoscopic simulation technique we used for HS 
particles in @ to colloids with several different attrac- 
tive interactions and also found that for strong enough 
attractions (but still below the point at which permanent 
clusters form) a similar maximum can be observed fl2l |. 
These simulations were performed at just one interme- 
diate Pc = 2.5, where the Peclet number is defined as 
Pc = iP s a/ Dq, a is the particle hydrodynamic radius, and 
Dq is the equilibrium self-diffusion coefficient for a single 
colloid at infinite dilution. 

Here we investigate the effect of varying the Pe number 
at one fixed attractive interaction strength. We find that 
the maximum in v s ((j)) v.s. (j) is stable for Pe < 1, but 
starts to diminish at Pc = 2.5 and has completely disap- 
peared by Pc = 10. We also measure the colloid- colloid 
radial distribution function g(r), and find that deviations 
from the equilibrium form begin to emerge for Pe > 1 . 

The properties of the clusters formed in irreversible ag- 
gregation processes under sedimentation has been stud- 
ied by means of Brownian Dynamics simulations in the 
region of low Pe number (Pc < 1) in 0,0]. Although 
permanent clusters don't form at equilibrium in our sys- 
tem, the attractions do promote transient cluster for- 
mation. Thus, fragmentation combines with aggregation 
and sedimentation processes studied previously. In con- 
trast to the Brownian Dynamics simulations, we also in- 
clude the effects of hydrodynamic interactions. We find 
that increasing the Pe number reduces the number of 
larger clusters, and that for a given cluster size, the num- 
ber of fragmentation events grows with increasing Pe, 
while the number of aggregation events decreases with 
increasing Pc. 

We proceed as follows: In the methods section we in- 
troduce our model potentials and describe the stochastic 



rotation dynamics computer simulation technique we em- 
ploy. In the results section we then present and discuss 
the effect of Pe number on v s (<p), g(r) and on cluster dis- 
tributions, and then summarize our findings in the final 
section. 



II. MODEL AND METHODS 

We employ a hybrid simulation method that treats the 
colloids with molecular dynamics (MD) and includes a 
background fluid that is by calculated with the Stochas- 
tic Rotation Dynamics (SRD) approach [lj|. Such a 
hybrid technique was first employed by Malevanets and 
Kapral [l6[ , and we adapted it to study steady state sed- 
imentation of HS particles in Q, and particles with at- 
tractive interactions in [l2[. Moreover, the method has 
recently been shown to quantitatively describe colloidal 
sedimentation experiments, including complex nonlinear 
effects such as Rayleigh Taylor instabilities (l7j. SRD 
itself (or Multi Particle Collision Dynamics as it is also 
called) has been fruitfully app lied to a large number of 
other problems as well [181 fl3 ] . These previous successes 
give us confidence to apply this method to the problem 
of Pe number effects on sedimentation. 

Within SRD, the fluid is not treated as a continuous 
solvent, but instead it is modelled by point-like particles. 
Thus the system is composed of N s solvent particles of 
mass m s , and N c colloids of mass M c . Collisions between 
solvent particles are efficiently coarse-grained in such a 
way that mass, momentum and energy are locally con- 
served, leading to the macroscopic Navicr-Stokes equa- 
tions (la. 118 - 1201 ] . At the same time, Brownian motion of 
colloids is achieved naturally through colloid-solvent col- 
lisions. As SRD also includes thermal noise, the method 
represents a powerful technique to simulate driven col- 
loidal suspensions embedded in a fluid. The choice of the 
simulation parameters must be done in a way that the 
most important dimensionlcss numbers fall in the cor- 
rect hydrodynamic regime. Care must also be taken to 
ensure that physically relevant time-scales are properly 
represented. We refer to Refs. |fj [2^ for further informa- 
tion about the choice of the SRD parameters. 

The interaction between solvent particles and colloids 
is represented by the repulsive WCA form (2~i| 



pV cs {r) 



r < 2 1 /V C 
r > 2 1 / 6 ct c 



(2) 

where r is the distance between the solvent particle and 
the center of the colloid, (3 = l/Zc^T, a cs is the colloid- 
solvent interaction range and a cs = 2.5. The use of radial 
colloid-solvent interaction implies that the fluid is not 
able to transfer angular momentum to the colloidal par- 
ticle. Consequently, it induces slip boundary conditions 
at the colloid surface. Although it is also possible to gen- 
erate algorithms with stick boundary conditions (22l424| , 
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we don't expect this to have any major qualitative affect 
on our results. 

The colloid-colloid interaction has been modelled by a 
classic DLVO [|[ potential, where the total interaction is 
the sum of three contributions, V cc = Vhs + VvcLW + Vdh- 
The first term is a repulsive hard sphere contribution 
given again by the WCA model but with a higher expo- 
nent to obtain a steeper hard-core repulsion 



PV cc (r) = 



r < 2 1 / 24 <r c 
r > 2V2V 



(3) 

where r now represents the distance between the centres 
of the colloids, e cc = 2.5 and a cc — 2A5a cs is the HS col- 
loidal diameter. The colloid-colloid diameter was slightly 
larger than 2a cs to avoid spurious depletion forces be- 
tween the colloids @, [20| • The second contribution is the 
repulsive Debye-Hiickel interaction between the electrical 
double layers, 

PV DH (r) = B— cxp[- K (r - a cc )} (4) 

r 



where B 



Z r is the effective 



charge of the colloid, e r and eo are the solution relative 
and vacuum dielectric constants respectively, and n is 
the reciprocal Debye screening length [25j]. Finally, the 
third term is the short-range London- van der Waals at- 
traction Q , given by. 



PV vd w{r) 



A 



H 



12k B T 



21n- 



(5) 

where Ah is the Hamakcr constant. In order to over- 
come the singularity of the van der Waals contribution 
at contact, we also introduced a cut-off distance, 5 (the 
so-called Stern layer) [26[ . 

We study two different limiting situations: hard-sphere 
colloids and attractive colloids with a short-range attrac- 
tion. In the former case, only the repulsive WCA poten- 
tial was used. In the latter case, the potential parameters 
were given by Z c = 8.3, k = 8.955/er cc , (3A H = 10.39 
and 5 = 0.048ct cc . The temperature was assumed to be 
T = 293K, and the relative dielectric constant of the sol- 
vent was given by the one for water, e r = 78.5. The re- 
sulting reduced second virial coefficient was B| = 1 .53. 
To calculate B^ , we divided the virial coefficient through 
that of a hard sphere with an effective diameter equal 
to the range of the short-range repulsive barrier of V cc . 
We note that there is some ambiguity in exactly how to 
obtain the effective diameter, as other recipes could be 
employed, but because of the steepness of the colloid- 
colloid interaction the differences are typically small, so 
we use this very simple one. 

The equations of motion are updated with a standard 
Molecular Dynamics algorithm for the colloid-colloid and 
the colloid-fluid interactions, and with the coarse-grained 
SRD collision step for the fluid-fluid interactions. The 



number density of solvent particles was p s = 40/cr^, the 
mass of the colloids M c = 168m s , and the stochastic ro- 
tation angle was fixed to a = ir/2. This particular choice 
has been shown to reproduce correctly the hydrodynamic 
behaviour in the low Reynolds number regime, as well as 
the thermal Brownian fluctuations and diffusion for col- 
loidal suspensions. The effective hydrodynamic radius 
obtained with these conditions is a « 0.8cr cs . We refer 
the reader to Refs. @, [2(| for further discussion of the 
technical details of our SRD parameters. 

The simulations were performed by placing N c = 
8 — 819 colloids inside an elongated box of sizes L x = 
L y = 16a cs and L z = 48<7 CS with periodic boundary con- 
ditions in the three spatial directions. The number of 
solvent particles was TV = 40Vf ree /a^ s ~ 3 — 5 x 10 5 , 
where Vf ree is the free volume left by the colloids. The 
sedimentation process is induced by applying a gravi- 
tational external field to the colloids in the z-direction. 
After an initial transient time, the system reaches steady- 
state conditions, where the average sedimentation veloc- 
ity v s is constant, the colloidal microstructure does not 
depend on time, and no drift is observed in other mea- 
sured obscrvablcs either. In our simulations, we var- 
ied Pe from 0.5 to 10 by increasing the gravitational 
field. Averages were collected over runs from 5000 Stokes 
times t s (for Pe = 0.5) to 50000t s (for Pc = 10), where 
ts = (T cs /Us is the time it takes an isolated sphere to 
settle one particle radius. One could also define a Brow- 
nian time tb = cr 2 s /Do, which is proportional to the 
time it takes the colloid to diffuse one radius. Then 
tb ~ TsPe, and averages are taken over periods rang- 
ing from 5000 to 10000 Brownian times, which should be 
ample for thermodynamic averaging. We work in a frame 
where the downward volume flux of colloids is compen- 
sated by an equivalent upward volume flux of fluid. In 
all cases the value of the Reynolds numbcri is kept small 
(Re = v s ajv < 0.4, where v is the kinematic viscosity), 
in order to ensure that the hydrodynamic behaviour of 
our system corresponds to the Stokesian regime typical 
for colloids. 



III. RESULTS FOR SEDIMENTATION 
VELOCITIES 

In Fig. [TK a) we compare v s (4>)/v® v.s. <j> for the HS and 
attractive colloid systems at four different Pe numbers 
(Pc = 0.5, 1, 2.5 and 10). The hydrodynamic packing 
fraction is given by <j> = ^irp c a 3 , where p c is the col- 
loidal number density. The hydrodynamic radius a is 
used because this sets the scale for the hydrodynamic in- 
teractions (one can also define a packing fraction based 
on colloid-colloid interactions, but this will typically be 
higher that the hydrodynamic one; the system is effec- 
tively non- additive [27j|). 

For hard spheres, as we increase the particle concen- 
tration, the hydrodynamic effect of the solvent leads to 
a rapid decrease of v s , as first predicted by Batchelor [BJ. 
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FIG. 1: (a) Normalized average sedimentation velocity 
(vs/Vg) as a function of the hydrodynamic packing fraction 
(4>) for HS and attractive colloids, calculated for Pe = 0.5, 1, 
2.5 and 10. (b) v s /v® for attractive particles at the low density 
regime. The dashed line represents the Batchelor prediction 
given by Eq. ([1} for the attractive system (B| = —1-53) in 
the limit Pe— > 0. The dotted line applies Batchelor's low den- 
sity limit equation, but employs the measured pair correlation 
function for Pe=10. 

As found in Q , v s (4>) /v® is virtually independent of Pe 
number, and is in good agreement with the theoretical 
prediction of Hayakawa and Ichiki [28| and with experi- 
mental data obtained for Pe — > oo. 

The picture changes dramatically when we consider 
our system with attractive interactions. On the one hand, 
for low Pe numbers, v s /v® exhibits non-monotonic be- 
haviour, with a maximum at intermediate 0, as predicted 
in [ll[ and observed by us in previous simulations [l2j]. 
According to Batchelor's pioneering calculations Q for 
the dilute limit, shown in Eq. (1), v s (<p) increases with <p 
if BjJ < —0.86, so our attractions are strong enough to in- 
duce an initial increase in v s (0) with </>. Indeed, as shown 
in Fig.QJb), at least for small Pe numbers (Pe < 2.5) the 
slope of the linear region at low <f>, i.e. the sedimentation 
coefficient, is consistent with the Batchelor theory that, 
through Eq. ((J), predicts that v s = v° « 1 + 2.3640, 
which agrees well with our simulation data. 

On the other hand, at high 0, the fluid must pass 
through a dense collection of colloids and so we expect 
v s ((j>) to eventually decrease with increasing <j). The 
crossover between these two regions explains the max- 



12 




FIG. 2: Different regions in the Pe-Bj diagram: v s (4>)/ v s 
either decreases monotonically with 0, or else shows a non- 
monotonic behavior (maximum). The diagram also has a 
third region for strong attraction where the colloidal sus- 
pension will phase separate and exhibit permanent clusters. 
E>2 < —1.5 is a good lower bound for where this behaviour 
will begin to set in [2!|. 



imum observed in Fig.QJa). 

For Pe < 1 the whole curve is independent of Pe, at 
least within the accuracy that we can measure. In con- 
trast to what was found for HS systems, increasing the Pe 
number above this regime does change the behaviour sig- 
nificantly. For example, there is a progressive reduction 
of the maximum. For Pe > 10, this reduction is strong 
enough to completely remove the maximum. The dilute 
regime also becomes affected. For Pe > 2.5 the sedimen- 
tation coefficient starts to noticeably decrease compared 
to the one predicted by the Batchelor. 

The critical value of the Peclet number that leads to 
this non-monotonic behaviour depends on the interac- 
tion potential, but for short-ranged interactions a max- 
imum can occur for sufficiently small Pe numbers if 
£?2 < —0.86. As Pe increases, we expect the maximum to 
disappear. Preliminary simulations at other B^s suggest 
that the weaker the attractions, the lower the critical Pe 
number at which the monotonic behaviour sets in. On 
the other hand, for strong enough attractions, perma- 
nent clusters will form and lead to different physics than 
what we are considering here. It is known that the crite- 
rion ^ 1-5 provides a remarkably accurate prediction 
of the potential strength at which a critical point sets 
in [2^| . Of course the critical point is only at one con- 
centration, so one could potentially reach a slightly more 
negative at lower concentrations as long as the fluid- 
fluid binodal or spinodal towards phase separation (and 
permanent cluster formation) is not crossed (note that 
the location of these lines will be affected by the sedi- 
mentation as well). One must also take care not to cross 
a fluid-solid binodal or other phase-lines. We illustrate 
the region of non-monotonic behaviour with a schematic 
diagram in Fig. [2] 
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FIG. 3: Radial distribution function g(r) as at different Pe 
number for HS (top) and attractive colloids (bottom). Note 
the differences in scale between the top and bottom graphs. 



IV. RESULTS FOR MICROSTRUCTURE 

In order to investigate further the possible causes for 
the change of v s (<f>)/v s with Pe, we also studied the mi- 
crostructure of the suspensions, including radial distri- 
bution functions and cluster distributions. 

Fig. [3] depicts the radial distribution function of the 
settling colloids at five different Pe numbers (from Pe = 
1 to 10) for HS and attractive colloids at <f> = 0.0233. 
For Pe < 1 the radial distribution functions for both HS 
and attractive colloids do not vary appreciably with Pe, 
and agree with the equilibrium distributions (at Pe = 
0). For larger Pe numbers we observe a small deviation 
for the HS correlation functions, including the emergence 
of a small but distinct hydrodynamically induced anti- 
correlation region around r « 1.8a cc for Pc = 10. Note 
that the scale of the changes is very small, so that the 
effect on v s ((j>) is expected to be very modest as well, 
most likely well within the error bars of our calculation 
of sedimentation coefficients we measure. 

In contrast to the case for HS particles, for attractive 
interactions the effect of Pe number on the pair correla- 
tions is much more pronounced. For Pe < 1, g C c(r) is 
indistinguishable from its Pc = form, which at these 
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FIG. 4: Average number of neighbours rineigh, given by 
Eq.©, around a colloid as a function of Pe for HS and at- 
tractive particles. 

low packing fractions can be accurately reproduced by 
the dilute limit expressions exp [— j3V cc {r)] As we 

increase the Pe number above 1, the main peak of g C c{f) 
begins to noticeably decrease in magnitude. At Pc=I0, 
the deviations are quite substantial. One measure of this 
hydrodynamic effect can be quantified by measuring the 
average number of nearest neighbours around a particle, 
which is given by the integral of g(r) from r = to the 
position of the first minimum(located at r w 1.5<r cc ): 

/*1.5<7 CC 

n n eig = 47rp c / g cc (r)r 2 dr. (6) 
Jo 

As shown in Fig. [H n ne i g drops from 2.39 to 2.02 as we 
increase Pe from 1 to 10. In other words, here the effect 
of increasing the relative strength of the convective over 
thermodynamic forces reduces the probability that two 
particles are in close proximity. 

In the p — > limit, instead of the equilibrium correla- 
tion function g(r) « exp(— (3V(r)), one could also use the 
measured non-equilibrium g(r) in the derivation leading 
to Eq. fl]). When we use the measured g(r) for Pe = 
10, we find that v s /v° s « 1 — 0.3260. The sedimenta- 
tion coefficient now has the opposite sign to that found 
when using the g(r) at Pe = 0, and in fact fits the data 
for P = 10 reasonably well, as can be seen in Fig. [TJb). 
This example illustrates the danger of using equilibrium 
correlation functions in non-equilibrium calculations. 

Next we focus on another measure of microstructure, 
namely cluster distributions. Two particles are consid- 
ered to be in a cluster if they are within a cutoff radius of 
r = 1 .06(T CC of each other, which places them well inside 
the attractive potential well. Other cutoffs could be used, 
but these don't qualitatively change the behaviour we ob- 
serve. The clusters are transient, that is they form and 
break up due to thermal and convective forces. Moreover, 
in steady-state, the average cluster distributions should 
not change with time (although fluctuations may occur 
caused by finite size effects). Using this criterion, we 
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FIG. 5: The plots show (a) the probability of finding a cluster 
of size i and (b) the average cluster life-time (normalized by 
the Brownian time tb) for attractive particles at several Pe 
numbers. In all cases the packing fraction of colloids is (j> — 
0.0233. 
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FIG. 6: Probability of aggregation and fragmentation events 
calculated for attractive particles at two Peclet numbers, Pe = 
1 and 10. The packing fraction of colloids is again cj> = 0.0233. 



measured the probability distribution of transient clus- 
ters P(i) and their average life-time Tj as a function of 
the cluster size, i. We used the same criterion to calculate 
cluster distributions for HS particles. 

For HS particles, both properties are only weakly af- 
fected by Pe number over the range we study. In con- 



trast, for attractive particles, Fig.[5ja) demonstrates that 
increasing Pe leads to a substantial decrease of the prob- 
ability of finding large clusters. At low packing frac- 
tions, this may help explain the decrease of the average 
sedimentation coefficient with increasing Pe number, as 
larger clusters are expected to sediment faster. 

Given that there are fewer large clusters with increas- 
ing Pe number, and that convective forces increase with 
Pe, one might also expect a decrease of the the cluster-life 
times as Pe grows. However, as can be seen in Fig. [5jb), 
for <fi = 0.0233, the cluster life-times are nearly indepen- 
dent of Pe number over the range of Pe numbers we inves- 
tigate. We note that this surprising effect seems limited 
to very low <p. At (f> ~ 0.05, we do observe a decrease in 
cluster life-times with increasing Pe. To explain why the 
average life-times of the clusters does not change at low </> 
it is important to think further about the processes that 
lead to a certain cluster size. These include both aggre- 
gation processes, whereby smaller clusters coagulate into 
larger ones, and fragmentation processes where a larger 
cluster breaks up into smaller ones. The larger the clus- 
ter, the more different formation and breakup pathways 
that are available. The life-time of a cluster can decrease 
both because the rate of aggregation increases or because 
the rate of fragmentation increases. 

To study the statistics of fragmentation and aggrega- 
tion, we tracked, over a suitably long averaging time, 
the number of times a given cluster of size i decreased 
in size due to fragmentation events, nf rag (i) as well as 
the number of times it increased in size due to aggre- 
gation events, n agg {i). These numbers increase linearly 
with time, so it is useful to define normalized quanti- 
ties P agg (i) = n agg (i)/n(i) and P frag (i) = n frag {i)/n(i), 
where n(i) is the total number of times that a cluster of 
size i emerges during the averaging time. In steady state, 
Pagg{i) + Pfrag{i) = 1, so that each can be viewed as a 
probability that a cluster of size i disappears by cither 
fragmentation or aggregation respectively. These proba- 
bilities are plotted in Fig.^for Pe=l and Pe=10. We see 
firstly that the probability of aggregation decreases with 
increasing cluster size, whereas the rate of fragmentation 
increases with increasing i, as one might expect because 
there are more ways a large cluster can break up. More- 
over, these probabilities change with Pe number. For 
higher Pe number the probability of fragmentation in- 
creases, as one might expect due to the increased rela- 
tive strength of shear forces over thermal forces. How- 
ever, the probability of aggregation also decreases due to 
the flow. The combination of these two rates apparently 
balance each other out, and so the average life-time of 
a cluster at <j) = 0.0233 appears not to change with in- 
creasing Pe number. At higher packing fractions these 
cluster life-times do change. By contrast, the increase 
in the probability of fragmentation, combined with the 
decrease in the probability of aggregation does lead to 
a decrease in the number of larger clusters. Clearly the 
effect of changing the balance of convective to thermal 
forces on the aggregation and fragmentation of clusters 
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is complex. A fuller investigation would not only track 
the cluster size distributions, but also the different path- 
ways of formation/breakup. A future publication will 
investigate these effects in more detail. 



V. DISCUSSION 

In summary, we have used a mcsoscopic simulation 
technique to study the effect of changing the balance be- 
tween thermal and convcctivc forces, quantified by the Pc 
number, in a scdimcnting suspension of colloidal parti- 
cles. Whereas for pure HS particles the effects of varying 
the Pe number between 0.1 and 10 are quite small @, for 
attractive particles there is a clear change in behaviour 
over this range of Pe numbers. Firstly, the maximum in 
the sedimentation velocity v s {4>) with <f>, first predicted 
in [ToL [Tl| , and observed by simulations in [l2j , disap- 
pears as the Pe number increases. Secondly, the radial 
distribution function g(r) changes noticeably for Pe > 1: 
increasing the Pe number means the average number of 
nearest neighbours decreases compared do the equilib- 
rium value. Thirdly, the cluster distributions change: 
the probability of observing larger clusters decreases with 
increasing Pe. Moreover, increasing the Pe number in- 
creases the probability that a transient cluster fragments, 
and decreases the probability that it aggregates with an- 
other particle or with other clusters. This work illus- 
trates the often complex crossover from equilibrium-like 
to purely non-equilibrium behaviour as the balance be- 
tween thermal and convective forces changes. 



A number of further questions are raised by this study. 
First of all, we only used one fairly short-ranged poten- 
tial form. It would be interesting to see how these re- 
sults depend on the shape of the potential. For exam- 
ple, a longer-ranged attractive potential with the same 
integrated strength, which through Eq. (TTJ) would gener- 
ate a similar sedimentation coefficient, would have a less 
deep dimcnsionless well depth Vmin/ksT. One might 
expect that finite Pe number effects set in roughly when 
Vmin/kBT ~ Pe, i.e. when the shear forces are strong 
enough to dislodge particles that cluster together [3l| . 
Thus a for a longer-ranged potential, finite Pe number 
effects should set in earlier than for a shorter-ranged po- 
tential. Similarly, it would be interesting to investigate 
what the effect of Pe number is on longer-ranged repul- 
sive potentials 32 1. 



Another interesting direction of investigation would be 
to study what happens at larger interaction strengths 
where permanent clusters can form. The shape of clus- 
ters [33J will most likely also depend on Pe number, and 
we anticipate a rich physics as a function of the interplay 
between the aggregation of clusters and the convective 
forces driven by their sedimentation. 
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